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ABSTRACT 


Supersonic  flow  past  two  oscillating  airfoils  is 
analyzed  in  the  thesis  using  an  analytical  elementary 
theory  valid  for  low  frequencies  of  oscillation.  The 
airfoils  may  have  arbitrary  stagger  angle.  This 
approach  generalizes  Sauer's  solution  for  a  single 
airfoil  oscillating  at  small  frequencies  in  an 
unbounded  supersonic  flow. 

It  is  shown  that  this  generalization  can  provide 
an  elementary  theory  for  supersonic  flow  past  two 
oscillating  airfoils.  This  aerodynamic  tool  will 
facilitate  the  calculation  of  pressure  distribution  and 
consequently  the  calculation  of  moment  coefficient. 
Torsional  flutter  boundaries  are  computed.  The  results 
for  the  pitch  damping  coefficient  are  the  same  when 
compared  with  previous  analysis.  For  arbitrary 
frequencies  a  linearized  method  of  characteristics  was 
outlined. 

The  elementary  theory  that  has  been  developed  in 
the  thesis  can  be  used  for  flutter  evaluation  of 
aircraft  carrying  external  stores.  The  result  of  the 
thesis  is  the  derivation  of  the  pitch  damping 
coefficient  which  is  necessary  to  predict  the  flutter 
conditions. 
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I.  INTRODUCTION 


Aeroelasticity  is  the  study  of  the  effect  of  aerodynamic  forces  on  elastic  bodies.  One 
of  the  interesting  problems  in  aeroelasticity  is  the  stability  of  a  structure  in  wind.  Since  for  a 
given  configuration  of  the  elastic  body,  the  aerodynamic  force  increases  rapidly  with  die  wind 
speed  while  the  elastic  stiffness  is  independent  of  the  wind,  there  may  exist  a  critical  wind 
speed  at  which  the  structure  becomes  unstable.  Such  instability  could  create  excessive 
deformation  and  failure  of  the  structure.  A  particular  problem  is  the  FLUTTER  of  structures, 
which  is  a  self-excited  vibration  phenomenon. 

Flutter  is  caused  by  the  interaction  of  all  three  forces:  aerodynamic,  elastic  and  inertial. 
Consider  as  an  example  a  cantilever  wing  mounted  in  a  wind  tunnel  (with  the  root  rigidly  built 
in).  Suppose  the  wing  is  deliberately  deflected  and  then  released.  At  sufficiently  low  wind 
speed,  the  oscillation  that  follows  this  disturbance  will  quickly  die  out.  However,  at  some 
higher  wind  speed  one  observes  a  steady  oscillation,  which  maintains  itself.  This  is  the  critical 
flutter  speed  for  the  wing.  At  higher  speeds,  the  oscillation  will  be  rapidly  divergent  causing 
structural  failure. 

When  large  external  bodies  such  as  engine  nacelles,  fuel  tanks,  or  electronic  warfare 
pods  are  added  to  the  wing  of  an  aircraft,  the  dynamic  characteristics  of  an  aircraft  will  be 
changed.  Over  the  last  several  decades,  there  has  been  considerable  interest  in  the  calculation 
of  the  unsteady  aerodynamic  forces  on  aircraft.  Many  of  the  methods  presently  in  use  for 
numerical  calculations  involve  application  of  a  two-dimensional  theory  to  solve  three- 
dimensional  flow  problems.  Most  of  these  methods  do  not  have  the  ability  to  model  an 
arbitrarily  shaped  aircraft.  However,  with  certain  assumptions,  it  is  possible  to  simplify  the 
flow  problems  to  adequately  account  for  the  effects  of  unsteady  aerodynamic  forces  on  aircraft 
components.  Unfortunately,  the  transonic  flight  regime  still  makes  it  very  difficult  to  determine 
the  unsteady  aerodynamic  forces  because  of  the  nonlinear  shock  and  flow  separation  effects, 
which  occur  in  this  flight  regime. 

A  good  example  for  the  requirement  of  the  wing / stores  flutter  analysis  is  the  fighter 

aircraft 
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Fighter  aircraft  are  commonly  designed  for  a  primary  mission,  such  as  air  superiority, 
which  may  require  few,  if  any,  wing-mounted  external  stores.  However,  for  increased 
effectiveness  and  versatility,  many  secondary  missions  evolve  which  are  necessary  and  require 
the  use  of  a  variety  of  external  stores.  Thus,  many  combinations  of  these  external  stores  must 
be  carried  at  various  stations  on  the  wings  to  achieve  the  complex,  multi-role  missions  required 
by  the  operational  commands.  Sometimes,  the  wing  stores  are  not  available  to  the  aircraft 
manufacturer  even  during  the  production  of  an  aircraft  and  in  that  case,  the  flutter  analysis 
should  be  conducted  using  data  not  available  from  the  original  design. 

A  typical  modem  tactical  fighter  can  carry  a  great  variety  of  wing  stores  and 
consequently  the  total  number  of  possible  aircraft/ store  configurations  is  huge.  Therefore,  on 
military  aircraft  the  effect  of  unsteady  aerodynamics  is  computed  for  both  the  clean  wings  (no 
stores)  and  the  wing  with  tip  missiles.  Using  these  two  aerodynamic  configurations,  the  flutter 
analysis  is  probably  done  on  between  300  and  400  selected  wing/ store  configurations.  The 
cost  for  the  above  effort  would  increase  by  several  orders  of  magnitude  if  wing/store 
aerodynamics  were  also  considered  for  each  of  these  configurations.  In  many  cases,  this  is  not 
necessary  since  the  store  aerodynamics  has  a  small  effect  on  the  flutter  speed.  However,  there 
are  cases  where  neglecting  the  store  aerodynamics  will  lead  to  the  overestimation  of  the  flutter 
speed.  Therefore,  the  determination  of  the  cost  effectiveness  of  such  analysis  is  a  critical  step. 

In  this  thesis,  we  first  review  the  most  important  papers,  which  have  been  published  in 
recent  years  to  analy2e  wing-store  flutter  problems.  It  is  evident  from  the  literature  survey  that 
most  approaches  are  based  on  modem  CFD  (Computational  Fluid  Dynamics)  approaches 
which  require  extensive  computations. 

However,  there  is  still  some  merit  in  simplified  analyses,  which  reveals  the  major 
physical  effects.  Limiting  our  approach  to  low  supersonic  flight  speed  we  propose  to  replace 
die  multiple  stores  underneath  the  wing  by  a  second  wing  and  to  apply  to  a  2-D  problem  the 
so  called  “  strip  theory  approach”  which  reduces  the  three  dimensional  flow  problem.  This 
simplification  has  the  advantage  of  accounting  for  the  interference  effects  between  two  lifting 
surfaces.  Moreover,  these  interference  effects  will  be  stronger  than  those  caused  by  the  actual 
stores.  The  insights  which  can  be  attained  with  such  an  analysis,  therefore  will  represent  an 
extreme  limiting  case,  but  it  will  nevertheless  reveal  important  physical  effects. 
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The  importance  of  such  interference  effects  is  well  known  in  transonic  and  supersonic 
wind  tunnel  testing  and  in  the  operation  of  transonic  and  supersonic  compressors.  Platzer 
(1973)  presented  an  analytical  theory  for  the  analysis  of  oscillatory  supersonic  wind  tunnel  and 
blade  interference  effects.  This  theory  is  based  on  the  assumption  that  the  actual  airfoil  or 
compressor  blade  can  be  replaced  by  a  flat  plate  and  that  it  oscillates  with  low  frequency.  This 
makes  it  possible  to  expand  the  solution  in  powers  of  the  reduced  frequency  and  to  retain  only 
the  zeroth  and  first  order  term.  Platzer  (1973)  showed  that  this  theory  provides  a  convenient 
analytical  way  to  estimate  the  pitch  damping  as  a  function  of  supersonic  Mach  number  and 
pitch-axis  (elastic  axis)  position. 

In  this  thesis,  we  modify  Platzer’s  solution  to  the  case  of  two  airfoils,  which  are  in 
close  proximity  to  each  other.  Two  cases  need  to  be  distinguished,  namely  the  one  with 
supersonic  leading  edge  locus  and  the  case  with  subsonic  leading  edge  locus.  These  two  cases 
and  the  basic  theory  will  be  explained  in  chapter  III.  The  extension  of  Platzer’s  theory  will  be 
given  in  chapter  IV.  This  is  followed  by  a  chapter  on  the  method  of  characteristic  approach, 
which  can  be  used  to  compare  with  the  analytical  results.  The  final  chapter  is  devoted  to  the 
flutter  analysis. 
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II.  LITERATURE  SURVEY 


EVALUATION  OF  METHODS  FOR  PREDICTION  AND  PREVENTION  OF 
WING/STORE  FLUTTER. 

Authors):  Pollock,  S.  J.  Sotomayer,  W.  A.;  Huttsell,  L.  J.;  Cooley,  D.  E 

Source:  Collect  Tech  Pap  AIAA  ASME  ASCE  AHS  Struct  Dyn.  Mater.  Conf.  22nd  AIAA 
Dyn.  Spec.  Conf.,  April  6-10  1981,  Atlanta,  GA,  New  York,  NY,  pp.  362-372 

In  response  to  the  need  to  reduce  costs  and  improve  safety  for  flutter  evaluation  of 
aircraft  carrying  external  stores,  the  Flight  Dynamics  Laboratory  (FDL)  has  sponsored  several 
efforts  in  the  technical  areas  of  unsteady  aerodynamics,  flutter  prediction,  and  active  flutter 
suppression.  This  paper  discusses  each  of  these  three  areas  as  they  relate  to  wing/ store  flutter 
and  presents  specific  examples  from  analyses  and  tests.  Steady  and  unsteady  pressure 
measurements  were  obtained  in  a  wind  tunnel  at  subsonic,  transonic  and  supersonic  speeds  on 
a  fighter  wing,  tip-mounted  launcher  and  store,  and  underwing  pylon  and  store.  Store  flutter 
calculations  were  performed  using  both  calculated  and  measured  data  to  determine  the 
influence  of  store  aerodynamics  on  the  flutter  characteristics.  To  improve  the  accuracy  and 
reduce  the  time  and  costs  of  flutter  evaluations  on  the  many  store  configurations  carried  by  Air 
Force  fighters,  the  Flight  Dynamics  Laboratory  (FDL)  has  sponsored  several  programs  in  the 
technical  areas  of  unsteady  aerodynamics  and  flutter  predictions.  Also,  the  FDL  has  been 
sponsoring  several  efforts  to  explore  the  potential  of  active  flutter  suppression  systems  using 
feedback  control  techniques  to  provide  the  required  stability  and  to  avoid  speed  placards. 

This  paper  reports  on  some  of  the  FDL  research  related  to  wing/store  flutter 
prediction  and  prevention.  The  research  includes  an  unsteady  aerodynamic  measurement 
program  for  a  representative  fighter  wing,  with  and  without  tip  missile  and  underwing  store 
with  test  data  covering  the  Mach  number  range  0.6  to  1.35.  Flutter  analyses  based  on  an  FDL 
computer  program  specifically  for  use  on  aircraft  with  external  stores  are  also  described. 
Flutter  trends  using  this  computer  program  are  presented  for  the  wing  with  and  without  stores 
based  on  sectional  force  coefficients  from  wind  tunnel  measurements  and  from  theoretical 
calculations.  A  brief  description  is  given  of  some  FDL  programs  in  active  flutter  suppression. 
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and  typical  results  are  presented  for  wings  with  stores,  which  indicate  significant  potential  for 
improvement  in  flutter  speeds. 

Flutter  trends  for  three  different  fighter  wing/  store  configurations  were  calculated 
using  the  FACES  flutter  analysis  procedure.  Use  was  made  of  available  modal  vibration  data 
and  measured  aerodynamic  data.  The  predicted  Mach  number  trends  gave  minimum  flutter 
speeds  at  transonic  and  low  supersonic  speeds  as  would  be  expected  from  the  trend  of 
measured  center  of  pressure  and  lift-curve  slope  data.  Although  the  aerodynamics  on  the  tip 
launcher  had  a  somewhat  detrimental  effect  on  flutter,  the  aerodynamics  on  the  tip  store  had  a 
much  larger  detrimental  effect.  For  the  underwing  store,  the  aerodynamic  effects  were 
beneficial  for  flutter  since  the  pylon  had  the  effect  of  decreasing  the  aerodynamic  loading  on 
the  outer  portion  of  the  wing. 


WIND  TUNNEL  TESTS  ON  A  FIGHTER  AIRCRAFT  WING/STORE  FLUTTER 

SUPPRESSION  SYSTEM 


Authors):  Turner,  M.  R. 

Source:  AGARD,  Rep.  on  a  Coop.  Program  on  Active  Flutter  Suppression,  Paper  Presented  at 
the  Structural  and  Material  Panel  Meet.,  50th  April  1980,  Athens,  Greece  AGARD. 

A  special  system  designed  using  analytical  data  was  tested  on  a  YF- 17  aircraft  model  in 
the  NASA  Langley  16'  wind  tunnel  and  succeeded  in  meeting  the  requirement  to  increase  the 
flutter  dynamic  pressure  by  70%  at  Mach  equals  0.  8.  The  system  was  designed  using  a  novel 
procedure,  which  provides  these  stability  margins,  uses  minimum  control  surface  movement  in 
turbulence,  and  can  be  designed  using  either  analytical  or  empirical  data.  Two  wing  tip 
accelerometers  and  a  leading  edge  control  surface  were  used.  Empirical  open  loop  transfer 
functions  obtained  during  die  test  showed  that  the  analytical  data  overestimated  the  response 
of  the  flutter  mode  to  leading  edge  control  surface  excitation. 
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AERODYNAMICS  USING  TIME-LINEARIZATION  TRANSONIC  FLUTTER 

ANALYSIS 


Authors):  Wong,  Y.S.,  Lee,  B.H.K.;  Murty,  H.S. 

Source:  Journal  of  Aircraft,  v.  30,  January-February  1993,  pp.  144-145 

A  survey  of  the  progress  made  in  the  development  of  numerical  simulation  techniques 
for  unsteady  transonic  flow  calculations  are  presented.  Computational  methods  in  three- 
dimensional  unsteady  transonic  flows  concentrate  mainly  on  the  transonic  small  disturbance 
equation  and  time-lineraization  approach.  An  algorithm  is  introduced  for  solving  flutter 
occurrence. 
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WING  FLUTTER  BOUNDARY  PREDICTION  USING  UNSTEADY  EULER 

AERODYNAMIC  METHOD 


Authors):  Lee-Rausch,  Elizabeth  M.;  Batina,  John  T. 

Source:  Journal  of  Aircraft,  v.  32, 2  March- April  1995,  AIAA,  Washington,  DC,  USA,  pp.  416- 
422 


Modifications  to  an  existing  three-dimensional,  implicit,  Euler/Reynolds-averaged 
Navier-Stokes  code  for  the  aeroelastic  analysis  of  wings  are  described.  These  modifications 
include  the  incorporation  of  a  deforming  mesh  algprithm  and  the  addition  of  the  structural 
equations  of  motion  for  their  simultaneous  time-integration  with  the  governing  flow  equations. 
This  article  gives  a  brief  description  of  these  modifications  and  presents  unsteady  calculations 
that  check  the  modifications  to  the  code.  Euler  flutter  results  for  an  isolated  45-deg  swept- 
back  wing  are  compared  with  experimental  data  for  seven  freestream  Mach  numbers  that 
define  the  flutter  boundary  over  a  range  of  Mach  number  from  0.499  to  1.14.  These 
comparisons  show  good  agreement  in  flutter  characteristics  for  freestream  Mach  numbers 
below  unity.  For  freestream  Mach  numbers  above  unity,  the  computed  aeroelastic  results 
predict  a  premature  rise  in  the  flutter  boundary  as  compared  with  the  experimental  boundary. 
Steady  and  unsteady  contours  of  surface  Mach  number  and  pressure  are  included  to  illustrate 
the  basic  flow  characteristics  of  the  time-marching  flutter  calculations  and  to  aid  in  identifying 
possible  causes  for  the  increase  in  the  computational  flutter  boundary. 
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APPLICATION  OF  TRANSONIC  SMALL  DISTURBANCE  THEORY  TO  THE 

ACTIVE  FLEXIBLE  WING  MODEL 


Authors):  Silva,  Walter  A.;  Bennett,  Robert  M 

Source:  Journal  of  Aircraft,  v.  32,  January-February  1995,  ALAA,  Washington,  DC,  USA,  pp. 
16-22 


A  code,  developed  at  the  NASA  Langjley  Research  Center,  is  applied  to  the  active 
flexible  wing  wind-tunnel  model  for  prediction  of  transonic  aeroelastic  behavior.  A  semispan 
computational  model  is  used  for  evaluation  of  symmetric  motions,  and  a  full-span  model  is 
used  for  evaluation  of  antisymmetric  motions.  Static  aeroelastic  solutions  using  the 
computational  aeroelastidty  program-transonic  small  disturbance  are  computed.  Dynamic 
(flutter)  analyses  are  then  performed  as  perturbations  about  the  static  aeroelastic  deformations 
and  presented  as  flutter  boundaries  in  terms  of  Mach  number  and  dynamic  pressure.  Flutter 
boundaries  that  take  into  account  modal  refinements,  vorticity  and  entropy  corrections, 
antisymmetric  motions,  and  sensitivity  to  the  modeling  of  the  wingtip  ballast  stores  are  also 
presented  and  compared  with  experimental  flutter  results. 


10 


WING-STORE  FLUTTER  ANALYSIS  OF  AN  AIRFOIL  IN  INCOMPRESSIBLE 

FLOW 


Authors):  Yang,  Zhi-Chun  Zhao,  Ling-Cheng 

Source:  Journal  of  Aircraft,  v.  26,  6  June  1989,  pp.  583-587 

The  flutter  of  two-dimensional  airfoil  with  external  store  is  analyzed  to  investigate  the 
effects  of  pylon  stiffness  on  flutter  speed.  Among  the  40  configurations  studied,  five  were 
tested  in  the  wind  tunnel  to  verify  the  analytical  results.  The  variations  of  wing-store  flutter 
speed  with  the  pylon  stiffness  can  be  divided  into  three  types.  The  curves  of  the  normal  and 
flutter  frequencies  vs  pylon  stiffness  have  the  same  pattern.  They  can  be  sketched 
approximately  by  the  aid  of  the  normal  frequencies  of  the  two  degenerated  two-degree-of- 
freedom  systems,  i.e.,  and  the  freely  hinged  and  rigidly  connected  store  cases.  A  limiting  flutter 
speed  for  very  small  pylon  stiffness  is  deduced,  which  is  useful  to  identify  which  type  of  flutter 
the  configuration  studied  belongs  to. 
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STUDY  OF  THE  EFFECT  OF  STORE  AERODYNAMICS  ON  WING/STORE 

FLUTTER. 


Authors):  Turner,  C.  D. 

Source:  Collect  Tech.  Pap.  ALAA  ASME  ASCE  AHS  Struct  Struct  Dyn.  Mater.  Conf.  22nd 
AIAA  Dyn.  Spec.  Conf.,  April  6-10, 1981,  and  1981,  Atlanta,  GA,  ALAA  (CP  812),  New  York, 
NY,  pp.  343-351,  Paper:  81-0604 

This  study  represents  the  first  systematic  analytical  study  of  the  effect  of  store 
aerodynamics  on  wing/store  flutter.  A  large  number  of  wing/store  single  carnage 
configurations  and  parameters  were  included  in  the  study;  multivariate  analysis  techniques 
were  used  for  the  first  time  to  analyze  wing/store  configurations,  modal  data,  and  flutter 
results.  The  results  of  the  multivariate  analysis  indicate  that  it  may  not  be  possible  to  develop 
general  guidelines,  but  it  is  possible  to  develop  specific  guidelines  for  use  with  a  particular 
aircraft.  This  study  was  the  first  attempt  to  do  a  systematic  analytical  study  of  the  effect  of 
store  aerodynamics  on  wing/ store  flutter.  To  determine  this  effect  flutter  analyses  were  done 
on  four  aircraft  with  single  carriage  of  three  basic  store  types.  In  all  308  configurations  were 
analyzed  with  and  without  store  aerodynamics  (tip  missiles,  tip  tanks,  and  underwing  stores). 
The  effects  of  stores  have  been  analyzed.  The  results  of  the  factor  analysis  indicate  that  it  may 
not  be  possible  to  develop  general  guidelines,  but  it  is  possible  to  develop  specific  guidelines 
for  use  with  a  particular  aircraft  The  conclusions  of  the  study,  as  far  as  it  concerns  the  effect 
of  store  aerodynamics  on  wing-store  flutter,  show  that  60%  of  the  tested  configurations 
require  a  change  in  the  flutter  speed  no  more  than  7%.  Therefore,  there  is  no  need  to  evaluate 
the  requirements  for  aerodynamic  modeling  of  the  store.  However,  for  the  rest  of  the 
configurations  there  is  such  a  need.  Moreover,  75%  of  the  last  reveal  that  nonrealistic  flutter 
results  are  obtained  without  taking  into  consideration  the  store  aerodynamic.  Of  course,  the 
results  of  this  study  represent  the  configurations  that  were  used  to  generate  the  data  base. 
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AERODYNAMIC  MODELING  OF  AN  OSCILLATING  WING  WITH  EXTERNAL 

STORES 

Authors):,  Sotomayer  W.  A.,  Dusto  A.  R,  Epton  M.A.,  Johnson  F.  T. 

Source:  AIAA  ,  New  York,  NY,  pp.  243-252,  Paper:  81-0609 


An  analysis  of  the  steady  and  unsteady  aerodynamic  forces  acting  on  a  fighter  aircraft 
wing  with  stores  has  been  done.  Computations  were  performed  with  paneling  methods 
capable  of  presenting  arbitrary  aircraft  configurations  in  subsonic  and  supersonic  flow. 
Interference  effects  from  a  tip  store  and  an  underwing  pylon  store  in  varying  stages  of 
completion  were  also  analyzed.  Detailed  comparisons  between  experimental  data  and 
numerical  computations  are  also  made. 

In  order  to  calculate  the  pressures,  forces  and  moments  on  an  aircraft  a  numerical 
approach  is  proposed.  Numerical  computations  are  based  on  solving  an  integral  equation 
formulation  of  the  flow  problem  being  considered.  Steady  state  calculations  were  done  with 
two  numerical  methods.  The  first  of  these  is  a  pilot  code  developed  by  Johns.  In  this  method, 
distributions  of  linearly  varying  sources  and  quadratically  varying  doublets  are  used  to 
represent  the  aerodynamic  surfaces.  Wakes  are  represented  by  doublets  with  constant 
streamwise  strength  and  linear  variation  in  the  spanwise  direction.  A  method  developed  by 
Woodward  was  also  employed  in  making  numerical  calculations.  This  method  utilizes 
distributions  of  sources  and  vortices  to  represent  an  aerodynamic  surface  and  the  wakes  shed 
from  various  components. 

For  unsteady  flow,  each  of  the  methods  makes  use  of  sources  and  doublets  as  the  basic 
aerodynamic  singularities.  In  a  method  developed  by  Johnson  for  unsteady  subsonic  flow, 
distributions  of  sources  and  doublets  are  used  to  represent  the  aerodynamic  surfaces.  In  the 
wake,  doublets  of  fixed  strength  and  position  are  used  to  represent  the  effects  of  unsteady 
wake  motion.  Calculations  were  also  done  with  a  special  doublet  lattice  method.  In  this 
method,  a  wing  is  represented  as  a  sheet  of  doublets,  components  such  as  the  fuselage  or  a 
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store  or  nacelle  can  be  represented  with  an  axisymmetric  distribution  of  sources  with  doublets 
on  the  surface. 

Approximations  are  imposed  on  these  distributions  and  on  the  kernel  functions  used  to 
represent  their  influences  on  the  flow.  These  approximations  are  briefly  described  along  with 
the  flow  boundary  conditions  imposed  on  the  aerodynamic  surfaces  and  wakes. 

It  is  interesting  that  the  small  disturbance  partial  differential  equation  is  being  used  for 
unsteady  in  viscid  incompressible  flow  along  with  the  linearized  form  expressed  in  terms  of  the 
velocity  potential  (p.  In  addition,  simple  harmonic  motion  is  assumed  as  in  the  thesis. 

(p=Re[(p*el5]  where  (p*  is  a  complex  quantity. 

Moreover  using  the  Helmoltz  equation  a  more  simple  relationship  is  established  for  <p*  with 
the  form  of  an  integral  equation.  An  approximate  solution  of  the  last  equation  permits  an 
evaluation  of  the  entire  flowfield.  Numerical  solutions  of  the  last  equation  are  achieved 
through  an  aerodynamic  influence  coefficient  method.  The  boundary  conditions  are 
established  taking  into  consideration  the  configuration  of  the  F-5  aircraft,  a  relatively  old 
aircraft.  The  aircraft  is  comprised  of  the  following  components: 

1)  wing 

2)  underwing  pylon 

3)  missile  rack 

4)  missile  body. 

As  in  the  thesis,  the  above  surfaces  are  approximated  as  infinitesimally  thin  surfaces  so 
the  thin  wing  theory  applies  except  the  for  missile  body.  The  boundary  conditions  are 
expressed  in  terms  of  Taylor  series  expansion  about  reference  surfaces  taking  into 
consideration  the  thickness. 

Initially  boundary  conditions  for  steady  state  flow  are  made  use  of  and  unsteady 
boundary  conditions  are  taken  as  a  perturbation  about  a  steady  state  condition.  In  addition, 
variations  in  Mach  number,  frequency  of  oscillations,  and  interference  effects  arising  from 
component  build  up  are  presented. 
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For  the  wing  missile,  the  boundary  conditions  are  of  mixed  type  (Dirichlet  and 
Newman).  On  the  surface  of  the  missile,  the  boundary  conditions  are  similar  to  the  other  parts 
of  the  aircraft. 

Using  a  Woodward  method,  a  lifting  surface  is  divided  into  numbers  of  aerodynamic 
panels  number,  containing  distributions  of  sources  and  vortices.  A  source  distribution 
represents  a  fuselage  [or  a  pod  or  external  store].  Wing  thickness  is  represented  by  a  lineady 
varying  source  distribution  in  which  the  strength  is  equated  to  the  chord  wise  slope  of  the  wing 
thickness.  Camber,  twist,  and  lifting  effects  are  represented  by  a  linearly  varying  vortex 
distribution  where  the  strength  is  determined  to  satisfy  tangential  flow  at  panel  control  points. 
An  iterative  procedure  is  employed  in  solving  the  boundary  value  problem.  For  analysis 
problems,  the  surface  slope  is  described,  and  singularity  strengths  are  determined  by  inverting 
the  matrix  of  aerodynamic  influence  coefficients.  With  the  strengths  of  the  aerodynamic 
singularities  known,  the  u,  v,  and  w  velocity  components  at  a  given  point  may  then  be 
determined.  Pressures,  forces,  and  moments  are  calculated  by  numerical  integration. 

In  the  doublet  lattice  method,  the  aerodynamic  surface  is  subdivided  into  a  series  of 
infinitesimally  thin  panels.  Along  the  quarter  chord  line  of  each  of  these  panels  is  contained  a. 
distribution  of  acceleration  potential  doublets.  The  strength  of  the  distribution  is  constant  but 
is  not  known.  Specification  of  the  normal  velocity  at  a  set  of  points  on  the  surface  determines 
the  loading  on  each  element  Locating  the  lifting  elements  at  the  quarter  chord  and  the 
collocation  point  of  the  three-quarter  chord  of  the  midspan  of  each  panel  usually  results  in 
reasonable  success.  An  approximation  for  the  lifting  pressure  coefficients  in  terms  of  the 
induced  normal  velocities  may  then  be  found.  lift,  moment,  roll  and  generalized  forces  may 
then  be  calculated. 

Recently,  a  series  of  wind  tunnel  tests  were  conducted  by  NLR  of  the  Netherlands  with 
sponsorship  from  the  Air  Force.  These  tests  involved  measurements  of  steady  and  unsteady 
pressures  and  forces  on  a  model  of  the  F-5  wing  with  external  stores.  Several  external  store 
arrangements,  which  are  described  and  shown  in  this  paper,  are  compared  with  numerical 
calculations.  Experimental  results  and  numerical  calculations  show  that 

1)  For  subsonic  flow,  the  interference  due  to  the  tip  store  is  much  greater  on  the 
outer  portion  of  the  wing  than  it  is  on  the  inner  portions. 


15 


2)  The  interference  effect  due  to  the  pylon  store  is  significant  on  the  lower  surface  of 
the  wing  for  subsonic  flow;  there  are  also  noticeable  effects  on  the  upper  surface  of  the  wing  at 
subsonic  speeds. 

3)  For  supersonic  speeds,  the  zone  of  influence  remains  localized  within  the  Mach 

cone. 

In  general,  numerical  calculations  and  experimental  results  agreed  reasonably  well  with 
each  other. 
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III.  PROBLEM  FORMULATION 


The  interference  between  the  airfoil  (wing  airfoil)  and  the  wing  stores  is  examined  by 
making  certain  assumptions  that  simplify  the  problem.  The  wing  airfoil  is  replaced  by  an 
oscillating  flat  plate  airfoil  and  the  wing  stores  by  another  oscillating  flat  plate  in  distance  d 
from  the  wing  airfoiL  Two  cases  need  to  be  distinguished,  namely: 


1)  Supersonic  flow  past  two  oscillating  airfoils  with  supersonic  leading  edge  locus 
shown  in  Fig.  1  (above) 

2)  Supersonic  flow  past  two  oscillating  airfoils  with  subsonic  leading  edge  locus 
shown  in  Fig.  2 
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The  differences  are  obvious  and  will  be  further  analyzed  during  the  examination  of 
each  case.  In  both  cases  multiple  reflections  occur.  As  explained  in  the  preface  the  interference 
between  the  two  airfoils  under  certain  conditions  could  result  in  unsteady  flow  phenomena  like 
flutter. 

Therefore,  the  study  of  the  oscillating  airfoils  is  important  in  predicting  the  flutter 
characteristics  or  dynamic  response  of  the  wing  airfoil  when  different  stores  are  attached.  For 
the  case  shown  in  Fig.  1: 

tan/?  <  cot  a  (III-la) 

Where  the  case  shown  in  Fig.  2  requires: 

tan  /?  >  cota  (Hl-lb) 

The  flow  is  assumed  a  non- viscous  compressible  two-dimensional  flow  of  a  perfect  gas 
governed  by  the  continuity  equation,  the  Euler,  and  the  energy  equation.  The  continuity 
equation  is: 


Dp  f  du  dv ' 


(ni-2) 


The  Euler  Equations: 

Du  ^  1  dp  _  q 
Dt  p  dx 

5r+I^= o 

Dt  pdy 

The  Energy  Equation: 

^  =  0 
Dt 


(III-3) 


(01-4) 


(HI-5) 
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Because  of  the  assumption  of  small  amplitude  oscillations  of  the  airfoils,  all  flow  quantities  are 
considered  small  perturbations  linearly  superimposed  on  free  stream  quantities.  The  velocities 
can  therefore  be  written  as: 

u  =um+u  (III-6) 

c  =  ca>+c  (ni-7) 

v  =  v  (in-8) 

In  addition,  the  pressure  and  density  perturbations  are  linearly  superimposed  on  freestream 
quantities.  Therefore: 


Ap-p-p^  (III-9) 

Ap  -p-px  (III-10) 


The  local  velocity  of  sound  is  given  by: 


c  = 


^  dp^2 

\dPj 


S=ct 


Furthermore: 

p 

-£—  -  const 
Pr 

Taking  the  total  differential: 

.ydp-,p~dp  =  0 


Or: 


(III-ll) 


(III-12) 


(III-13) 
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(ni-i4) 


dp_ 

dp 


So: 


(III-15) 


From  AA  4318  class  notes,  the  surface  boundary  conditions  will  be  calculated  Requirement  is 
that  no  flow  will  penetrate  the  solid  surface  of  the  airfoil  resulting  in  the  flow  tangency 
condition: 


(y-Vs)-n  =  0  (HI-16) 

where  V,VS  are  the  velocities  of  the  free  stream  and  the  surface  over  which  the  fluid  flowsJf 
the  body  is  described  by  the  equation: 

F(x,y,t)  =  0  (III-17) 


then  the  normal  to  the  surface  at  any  point  is  given  by  the  following  relation: 


-  VF 

n  =  -. - r 

|VF| 


(III-18) 


Taking  into  consideration  that  the  substantial  derivative  of  a  surface  particle  should  be  zero: 


dF_ 

dt 


+  Fi-VF  =  0 


one  obtains  from  (11-16): 


V-VF  =  VS-VF 


(III-19) 


(ni-20) 


and  therefore: 
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(III-21) 


DF 

Dt 


dF_ 

dt 


+  V-WF  =  0 


Or  for  two  dimensional  flow: 


DF  dF  dF  dF 

- =  —  +  u —  +  v —  = 

Dt  dt  dx  dy 


(III-22) 


Particles  in  contact  with  the  surface  must  have  the  same  normal  velocities  as  the  surface.  This 
means  that  the  rate  of  change  of  F  is  zero  hence: 

—  =  0  (in-23 ) 

Dt 


For  an  airfoil  the  equation  of  the  upper  surface  can  be  written  as  follow: 

Fu(x,y,t)  =  y-  yu  (x,  t)  =  0  (III-24) 


Where  yu  is  the  distance  from  the  chord  line  to  the  upper  surface.  The  equation  of  the  lower 
surface  can  be  written  in  a  similar  way: 


FL(x,y,t)  =  y-yL(x,t)  =  0 


Aty=yu: 


DFu  _ 
Dt 


_(M).^L  +  V  =  0 

dt  dx 


At  y=yL: 


DFj  dyr  -  dyi 

— A  = Ak_(w)._£2^  +  V  =  0 

Dt  dt  dx 

=  =  i 

dy  dy 


(III-25) 


(III-26) 


(III-27) 


(III-28) 
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V(*>yi>')=v(*’°v) 


(III-37) 


Therefore: 


(111-38) 


4ov)=^+K.^=r 


(III-39) 


For  simple  harmonic  motion  of  the  airfoil: 


yu  =  K(xYa 


(111-40) 


=  h(xY 


(III-41) 


Therefore,  the  linearized  flow  tangency  equations  can  be  written: 


v(x,0+,t)=  [iahu(x)  +  ux^]e“*  ,..y  =  0+ 


(m-42) 


v(x,0",/)=  [iayL  (x)  +  ,..y  =  O' 


(III-43) 


The  above  relationships  describe  the  normal  velocity  on  the  lower  airfoil  (store  airfoil). 
Similarly,  one  obtains  for  the  upper  airfoil  (wing  airfoil)  at  y=d: 


v{x,d',t)=  \icohL{x)  +  =  d~ 


(III-44) 


if  the  store  airfoil  executes  a  pitch  oscillation  about  x  —  x0  one  has: 


y  =  -0(x-xoy 


(01-45) 


v(x,0+,0  =  [-ux0  -  i(od{x  -  x^e'*  at  y  =  0+ 


(III-46) 
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and  similarly  for  the  wing  airfoil: 


y  =  -6(x  -B-  jt0>'W)  (PLAT) 

v(x,  cT  j)  -  [-8{u„  cos  S  -  a. sin  S(x  -  B  -  xQ }-  i  #{«„  sin  8  +  co.  cos  S(x  -  B  -  x0 )}]e,a" 


(III-48) 
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IV.  ELEMENTARY  THEORY 


The  case  of  supersonic  leading  edge  locus  will  be  examined  assuming  first  that  both 
airfoils  are  oscillating  slowly  in  a  supersonic  flow.  The  airfoils  execute  vibrations  of  identical 
modes  and  amplitude  but  with  different  phase  angle  8  between  them.  It  is  obvious  that  the 
flow  field  is  such  that  no  disturbance  can  propagate  upstream  of  the  airfoil  leading  edges. 
Consequently,  the  perturbation  potential  can  be  written  as: 

$(x,y,t)  =cp(x,ydj)eto  (IV-1) 

where  x,  y  are  non  dimensional  coordinates  and  t  is  non  dimensional  time.  The  differential 
equation  for  the  perturbation  potential  in  a  non-dimensional  form  then  is: 

cota-(pIX-(pyy+2ikM2cp -k  2M2(p=0  (TV-2) 

The  linearized  boundary  conditions  are: 

tpy(x,0)=v(x)  on  the  lower(store  airfoil)  at  y=0  (IV-3) 

where  v(x)  is  the  downwash  amplitude. 

.C 

(py(x,d)=-e  ‘  v(x-dtanP)  on  the  upper(wing)  airfoil  at  y=d  (TV-4) 

where  it  is  assumed  that  the  two  airfoils  oscillate  with  a  phase  angle  8. 

Assuming  oscillation  around  Xq  v(x)  becomes: 

v(x)=0[l +ik(x— Xq)]  and  without  loss  of  generality  we  can  set  0=1  (TV-5) 

For  small  oscillations  (low  frequency)  the  potential  cp  can  be  expanded  and  neglecting  the 
higher  order  terms  of  k  one  can  write: 

(p(x,yh)  =X(x,y)  +kvP  (x,y)  (EV-6) 
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Therefore  the  differential  equation  (II-2)  becomes: 


cot2a-X)n-Xyy=0  CTV-7) 

cota-'Fyy-'Fyy+  2iMX=0  (IV-8) 

A  general  solution  of  the  above  system  of  equations  has  been  given  by  Sauer  (1950).  It  reads: 

X(x,y)=g(z)  (TV-9) 

'F(x,y)=h(z)+iMyg(z)/cosa  where  z=x-ycota  (TV-10) 

'T(x,y)=  h(z) +iMy  g(z)/cosa  where  z=x+ ycota  (TV-11) 

h(z)  and  g(z)  are  arbitrary  functions  for  positive  arguments  of  z  and  zero  for  z<0,  it  is  obvious 
that  the  expressions  of  the  above  functions  should  be  consistent  with  the  last  assumption. 

Using  the  two  solutions  for  left  and  right  running  Mach  waves  in  the  supersonic  flow 
field  of  slowly  oscillating  airfoils  pressure  coefficients  can  be  obtained.  In  order  to  do  so  the 
flow  field  between  the  two  airfoils  has  to  be  divided  into  several  zones.  The  number  of  zones 
depends  upon  A.  It  is  obvious  that  for  A>1  where  A=dcota  there  is  one  zone  along  each 
airfoil,  for  0.5<A<1  there  are  two,  for  0.33<A<0.5  there  are  three  and  so  on.  In  the  first  zone 
I  there  is  no  interference  from  the  airfoil  wing  (upper  wing),  which  means  that  the  lower  airfoil 
(zone  I  only)  does  not  sense  the  upper  airfoiL  In  that  case,  the  solution  has  the  form: 

Cp(x,yjk)=g(z)+k  (h(z)—  iMyg(z)/ cosa} 

The  boundary  conditions  that  must  be  satisfied  at  y=0  are: 

XjP-cota-g'  (z)=-l 

Yy=-cotah'(z)-iMyg(z)/cosa=-i(x  -Xo  ) 

From  the  last  two  equations  it  is  obvious  that: 
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A=0.5 


SUPERSONIC  FLOW 
PAST  TWO 

OSCILLATING  AIRFOILS 


Fig.  3 
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g(z)=ztana 


h(z)=— iztana(x0+ztan2Ot/2) 

Since  we  are  assuming  a  slowly  oscillating  airfoil,  the  solution  to  the  PDE  has  the  form: 
(p(x,y^)=g(z)+k|h(z)-iMyg(z)/ cosa] 
where  z=x-ycota 

In  order  to  calculate  the  pressure  distribution  on  the  upper  surface  of  the  lower  airfoil  in  the 
zone  I  the  following  relation  should  be  used: 

cp=-2fik(p  +(pj 

Substituting  the  known  functions  cp  and  cpx 

<p=ztana+k{-iztana(x0+ztan2a/2)  -(lMyztana/  cosa)} 

(px=d/dx[ztana  +k{-iztana(x0+ztan2a/2)  — (iMyztana/cosa)}]= 

=tana-iktana[x0+ztan2a-(Myz/ cosa)]= 

=tana— ik[tana(x0+ztan2a)+(y/ sinacosa)(sina/  cosa)]= 

2  2 

=tana-ik[tana(x0+ztan  a)+(y/cos  a)]= 

=tana-ik[tana(x0+  (x-ycota)tan2a+ (y/ cos2a)] = 

=tana-ik[tana(x0+x  tan2a-ytana+(y/ cos2a)]= 

=tana— ik[tana(x0+  x  tan2a-ycota]= 

Therefore,  the  pressure  distribution  is  given  by  the  following  relationship: 
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cp=-2{tana-ik[tana(x0+xtan2a— ycota]  +ik[ztana+k{-iztana(x0+ztan  a/2)— iMyztana/ 
cosa}]} 

Taking  into  consideration  that  higher  order  terms  of  k  are  neglected  the  pressure  distribution  is 
given  by  the  following  relationship: 

cp=2tana{— l+ik[x„+x  (tan2a— l)+2ycota]} 

and  therefore  on  the  airfoil  for  y=0: 

cp=2tana{-l+ik[x0+x  (tan2a— 1)]} 

In  die  zone  I*  the  perturbation  potential  can  also  be  written  direcdy  from  Sauer’s  single  airfoil 
solution: 

<P(x,y,k)  =gi(zi)+k[h1(z1)+  iMyg(zi)/cosa] 

Where: 

Z!=x— A— B +ycota 

The  boundary  conditions  that  must  be  satisfied  at  y=d  are: 

Xy=cota-g'(z1)  =  -e‘5 

*8 

Ty  =cota-h/ (zj)— iMgj (Zi)/ cosa  +iMy  gi'(zj)=-iel  (x— Xo~B): 

After  integration  and  algebraic  manipulations: 

gi  (zi)=-zitana-e5 

h1(z1)=:iz1tana*e1^(xo+z1tan2ot/2+Md/ cosa) 
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where  again  the  condition  has  been  imposed  that  the  functions  gi(zi)  2nd  hi(zi)  are  zero  for 
negative  and  zero  arguments. 

Following  the  same  procedure  as  in  zone  I,  the  perturbation  potential  (p  will  be 
calculated  substituting  the  known  functions  gi(zi)  and  hi(zi): 

<p=— 2itanae*^+ik{zitana  e1^(x0+Zitan2ct/2+Md/cosa)— Myzjtanae  /cosa} 

(pI=— tanote'^+iktancx  e1^{(xo+Zitan2a+Md/ cosa)— My/cosa}  — 

=-tana-eiS+iktana-ei5  { (xo+(x-A-B+ycota)tan2a+Md/cosa)-My/ cosa}  = 

=-tana-ei5+iktana-ei5  {  (xo+  (x-B)tan2a-A  tan2a+ycota+d/sina-cosa) 

— y/sinacosa}  = 

=-tana-eiS+iktana-ei5  {  (xo+(x-B)tan2a+A-ycota} 


Because  d= A/cota ,  tana-1  /cosa-sina-cota 

This  gives  the  following  result  for  the  pressure  distributions  on  the  lower  surface  of  the  upper 
airfoil  (wing  airfoil)  in  zone  I*: 

cp=2tana-ei5-2iktana-ei5{(xo+(x-B)(tan2a-l)+2A-2ycota}= 
cp=2tana-ei5  {l-ik{  (xo+(x-B)(tan2a-l)+2A-2ycota} 

and  therefore  at  y=d: 

cp=2tana-e'^  { 1-ik  {  (xq+  (x-B)(tan2a-l)} 
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Assuming  that  0.5<A<1  in  addition  to  Tones  I,  I*  there  are  two  more  Tones  II,  II*.  In 
Tone  II  the  perturbation  potential  is  due  to  the  initial  waves  from  zones  I,  I*  plus  the  reflected 
wave  at  zone  ILSo  the  reflected  waves  must  be  determined  by  adding  the  functions  g>(z£  and 
h2(z2).  In  the  zone  II  Sauer’s  solution  for  a  single  airfoil  applies  and  following  the  same 

procedure  as  previously: 

<p(x,y,k)  =g(z)+k[h(z)-iM(yg(z)/ cosa  l+giCzO+kjh^zO+iMfygitzi)/ cosa]+ 
g2(22)+k[h2(Z2)-iMyg2(Z2  )/ cosa] 

The  variable  z  has  the  following  form: 
z2=x-2A+ycota 

This  ensures  that  z2<0  upstream  of  the  leading  edge  of  the  airfoil. 

The  boundary  conditions  that  must  be  satisfied  at  y=d  are: 

i8 

Xy=-cota*gr(z)+cota-gif(z1)+cota-g2f(z2) — e 

f 

^y=-cota-hi,(z)+cota-h1,(z1)+cota-h2,(z2)-iMlg(z)-gi(z1)-g2(z^]/cosa+iMdcotalgi(z1) 

t  t  jg 

+gi  (zi)+g2  (Z2 )]/ cosa=-ie  (x-xo-B) 

After  integration  and  algebraic  manipulations  similar  to  the  previous  calculations  for  zones  I. 
I*: 

g2(z2)=ztana 

h2(z2)  =-  iz2tana  (xq+  z2/2tan2a  +2Md/ cosa) 

Therefore,  the  perturbation  potential  is  given  by  the  following  relationship: 
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(t>I=2tana-tana-ei8-iktana(x0+z2tan2a+2Md/cosa)+iktana-e5(x0+z1tan2a+Md/cosa) 

i8 

-iktana(x„+  ztaa2a  )-ikMy[tana-tana+  tana-e  ]/cosa= 

=2tana-tana-ei5-iktana[(x0+z2tan2a+2Md/cosa)-ei6(x0+  Zitan2a  +Md/cosa)-(x0+  ztan  a 

)+Myel8  /cosa]= 

=2tana-tana-ei5-iktana[x0+(x-2A+ycota)tan2a+2Md/cosa-ei5(x0+(x-A-B+ycota)tan2a 

i5 

+Md/cosa)-(x0+(x-ycota)  tan2a)+Mye  /cosa]= 

=2tana-tana-ei6-iktana[x0+(x-2A)tan2a+2Atafla/sma-cosa-e,6(x0+(x-A-B+ycota)tan2a 

*8 

+Atana/sina*cosa)-(x0+x  tan2a-ye  /sina-cosa)]= 

'8 

=2tana— tana-ei8— iktana[x0+xt2n2a+2A-ei8(x0+(x-B+ycota)tan2a+A)+x0+xtan2a+yel 

/sina*cosa]= 

=2tana-tana*ei8-iktana[2x0+2xtan2a+2A-ei8(x0+(x-B+ycota)tan2a+A)+ye  /sincx-cosa]= 

=2tana-tana-ei5-iktana[(2-ei5  )(xo+xtan2a+A)+  ei5(Btafl2a+  ycota)]= 

=tana{2-ei8— ik[(2— ei8  )(xo+xtan2a+A)+ei5(Btan2a+  ycota)]} 

This  gives  the  following  result  for  the  pressure  distributions  on  the  lower  surface  of  the 
upper  airfoil  (wing  airfoil)  in  zone  II: 

Cp=-tana{2-ei8 +ik[ei8(x0+  (x-B)(tan2a-l)+2A-2ycota)- 
2(x0+x(tan2a-l)+2A)} 


and  for  y=d  we  have: 
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Cp=-tana{2— e1^+ik[e1^(x0+(x-B)(tan2a-l))-2(x0+x(tan2a-l)+2A)} 

In  zone  II*  the  perturbation  potential  is  due  to  the  initial  waves  from  zones  I,  I*  plus 
the  reflected  wave  at  zone  II,  H*.So  the  reflected  waves  must  be  determined  by  adding  the 
functions  §3(23)  and  113(23).  Therefore,  in  the  zone  II*  Sauer’s  solution  for  a  single  airfoil 

applies  and  following  the  same  procedure  as  previously: 

<p(x,y,k)  =g(z)+k[h(z)-iM(yg(z)/cosa]+gi(z1)+k[h1(z1)+iM(ygi(z1)/cosa]+ 

g3(z3)+k[h3(z3)-iMyg3(z3  )/cosa] 

The  variable  Z3  has  the  following  form: 

Z3=x— A— B-ycota 

The  boundary  conditions  that  must  be  satisfied  at  y=0  are: 

Xy  =-  cota-g'  (  z)  +cota-gi'  ( zi)  -cota-g3f  ( z3)=  -1 

vFy=-cota-h'(z)+cota-h1'(z1)-cota-h3,(z3)-iM|g(z)-g1(z1)+g3(z3)]/cosa=-  i(x-xo) 

After  integration  and  algebraic  manipulations  similar  to  the  previous  calculations  for  zones  I, 
I*: 


g3(z3)=-z3tana-e:8 

h3(z3)  =iz3tana  (xo+Z3tan2a/2+Md/ cosa)ei8 

Therefore,  flie  perturbation  potential  is  given  by  the  following  relationship: 

(pI=tana-2tana-e^+iktana{ei^(xo+z1tan2a+Md/cosa)+e^(x0+z3tan2a+Md/cosa)-(x0+ 

ztan2a  )— My/ cosa} 
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tS 

(px=tana-2tana*e^+iktana { e'^[xo+ (x-A-B+ycota) tan  cx+Md/coscx]e  [xo+ 

(x— A-B— ycota)tan2a  +Md/ cosoc]  — (xq+  (x— vcotcx)tan  a)— My/cosa}  — 

(pI=tana-2tana'ei5+ik{ei52tana[x0+(x-A-B)tan2a+Md/cosa]-tana(x0+(x-ycota)tan2a)- 

y/cos2a} 

(pI=tana-2tana-ei5+ik{ei52tana[xo+(x-B)tan2a-Atan2a+Adtana/sma-cosa]-tana(x0+ 

(x-ycota)tan2<x  )-y/ cos2a} = 

(pI=tana-2tana-ei5+ik{ei52tana[x0+(x-B)tan2a+A]-taoa(x0+(x-ycota)tan2a)- 

y/cos2a}= 

(pI=tana-2tana-ei5+ik{ei52tana[x0+(x-B)taQ2a+A]-tana(x0+xtan2a+  ycota)} 

Because:  y/cos2a+  ycotatan3a=— y 

In  order  to  calculate  the  pressure  distribution  on  the  upper  surface  of  the  lower  airfoil  in  the 
zone  II*  the  following  relation  should  be  used: 

cp  =-2[ikcp  +<px] 

Substituting  the  known  functions  <p  and  (px 

=2tana-4tana-ei5-2ik{ei52tana[xo+(x-B)tan2a+A]-tana(x0+xtan2a+  ycota)} 

-2ik{ztana-z,tana-  ei5-  z3tana-  ei5}+higher  order  of  k  terms  which  will  be  neglected 
Therefore: 

Cp=2tana-4tana-ei5-2ik{ei52tana[x0+(x-B)tan2a+A]-tana(x0+xtan2a+  ycota)} 
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O  .g 

— 2ik{  (x-ycota)tana-(  x-A-B+ycota)tana-  e1  — ( x-A-B-ycota)tana*  e  } 

After  several  relatively  simple  algebraic  manipulations: 

Cp=-2tana[l-2ei5-ik{2ei5[xo+(x-B)(tan2a-l)+2A]-(x0+x(tan2a-l)+2ycota)}] 

Since  tiie  pressures  have  been  calculated  for  the  flow  field  between  the  airfoils  a  similar 
approach  will  be  utilized  to  obtain  the  pressures  on  the  bottom  of  the  lower  airfoil.  In  that 
case,  the  solution  has  the  form: 

<P(x,y>k)  =  g  (2)  +k{h(z)  +iMyg(z)/cosa} 

The  boundary  conditions  that  must  be  satisfied  at  y=0  are: 


Xy  =cota-g'(z)=  -1 


*Py=  cota-h'(z)+iMg(z)/cosa=-i(x-Xo)  aty=0 
From  the  last  two  equations  it  is  obvious  that: 
g(z)=ztana 

h(z)=-iztana(x+ztan2a/2) 

The  solution  again  has  the  form: 
cp(x,y^)=g(z)+kpi(z)+iMyg(z)/cosa] 
where  z  =  x+ycota 

In  order  to  calculate  the  pressure  distribution  on  the  upper  surface  of  the  lower  airfoil  in  the 
zone  I  the  following  relation  should  be  used: 

cp  =-2[ik(p +<p  J 
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(p=-ztana  +k{-iztana(x+ztan2a/2)+(iMyztana/ cosa)} 

(p^d/dxj— ztana+k{— iztana(x+ztan2a/2)+(iMyztana/cosa)}]= 

=— tana—  iktana[x+ztan2a-(Myz/cosa)]= 

=-tana-ik[tana(x0+ztan2a)— (y/sinacosa)(sina/cosa)]:= 

2  2 

=tana  -  ik[tana(x  0+  ztan  a)+(y/ cos  a)]  = 

Therefore,  the  pressure  distribution  is  given  by  the  following  relationship: 

cp=2tana{  1— ik[tana(x0+  x  (tan2a-l)+2ycota} 

Using  the  results  for  zones  I,I*,II,II*  it  is  possible  to  calculate  the  pressure  jump  across 
the  bottom  airfoil  (wing  store  airfoil)  from  the  following  relationship: 

Acp(x,0)=cp(x,0+)-  cp(x,0”) 

Two  integration  intervals,  shown  in  Fig. 4,  can  be  distinguished  in  the  lower  airfoil  taking  into 
consideration  that  0.5<A<1: 

1)  0<x<A+B 

2)  A+B<x<2A 

In  the  first  integration  interval: 

Acp,(x,0)=cp(x,0+)-  cp(x,0> 

=2tana{  l-ik[tana(x0+x  (tan2a-l)+2ycota) 

+2tana{  l-ik[tana(x0+x  (tan2a-l)+2ycota}  and  y=0 
Acp](x,0)=4tana{  l-ik(xo+  x(tan2a-l)} 


In  the  second  integration  interval: 
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+2tana[l-2ei5+ik{2ei5[xo+(x-B)(tan2a-l)+2A]-(xo+x(tan2a-l)+2ycota)}] 
=+2tana[l-2ei5+ik{2ei5[xo+(x-B)(tan2a-l)+2A]-2(x0+x(tan2a-l))}]  since  y=0 

But  e^=cos8+isin5  and  the  last  expression  can  be  written: 

=+2tana[l-2(cos8+isin8)+ik{2  (cos8+isin8)[xo+(x-B)(tan2 a*-l)+2A]-2(x0+x(tan2a-l))}] 

This  expression  is  a  complex  equation  with  real  and  imaginary  part  In  order  to  calculate  the 
pitch-damping  coefficient  only  the  imaginary  part  of  Acp  is  needed.  Therefore  one  has. 

TmAcpI(x,0) =-4tana(x0+x(tan2(X— 1)) 
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ImAcpn(x,0)=+2tafla[-2sin8/k+2cos5[xo+(x-B)(tan2a-l)+2A]-2(xo+x(tan2a-l))] 

The  lift  and  moment  coefficient  can  be  calculated  from  the  following  relationship. 

Q^of  JAcpdx 
CM=ol  1Acp(x-x0)dx 
By  definition: 

CM=0[  C,^  +ik  Cj^e'] 

Qm0'=  0^1+ 

(1  /ik) JA+BImAcpI(x,0)  (x-x^cIx+a+bI  1ImAcpn(x,0)(x-xQ)dx= 

JA+B-4tana(x0+x(tan2a-l)(x-xfl)dx+A+BJ2tana[-2sin5/k+2cos5[xo+(x-B)(tan2a-l)+2A]-2( 

Xo+x(tan2a-l))]  (x-x0)dx 

The  last  two  integrals  will  be  calculated  separately.  The  first  one  is: 

2tana[-2xo2  (A+B)+  x0(A+B)(A+B)-2(tan2a-l)[x0(A+B)2/2-(A+B)3/3)]  (TV-12) 

The  second  one  can  be  calculated  using  Maple.  The  results  are  a  complicate  expression. 
Two  calculations  have  been  included  for  the  convenience  of  the  reader  since  in  both  integrals 

there  is  a  common  factor  2tana.  The  result  is: 

2tana{2x0+tan2a*B3A2B-cos8-tan2a-B3/3+2tan2a/3-sin5-A2/k-sin5-B2/k+sm8/k+2cos5/3 

+2x0sin8*A/k-2x02cos8-B-2x02cos8-A+4x0cos8-A+2x0cos8-B+2x0sm8*B/k+2x02A+2x02B 

+2cos8*tan2aA3/ 3+2cos8-AB2-4x0AB-2sin8-AB/k-2/3+4cos8- A3/3+2B3/3+2AB2+2B  A2 

-2tanza-AB2-2tan2a-AB2-2cos8-tan2a/3— 2x0cos8-x0tan2a-cos8<B-2cos8-A-2x0sin8/k+2x02 
sin8/k+2x02cos8-2x02+x0tan2a-A2+x0cos8-tan2a+x0tan2a-B2+cos8-Btan2a-2x0cos8-Btan2a+ 
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3cosS  A.2B-2tan2aA.3/ 3-2tan2a-B3/  3+cos8-B3/3+2A3/3-2x0B2+2x0tan2a-AB 
^x^^XoCosS-BA+XoCOsS-tan^-B^XoCosS-A^XoCosS-tan^-A2}  (IV-13) 

Therefore: 

CMe'=2tana[-2x02  (A+B)+  x0(A+B)(A+B)-2(tan2a-l)[x0(A+B)2/2-(A+B)3/3)]+ 
2tana{2x0+tan2a-B3A2B-cos8’tan2a-B3/ 3+2tan2a/3— sin8-A2,/k— sin8-B2/k+sin8/k+2cos8/3 
+2x0sin8-A/k— 2xo2cos8-B— 2x02cos8-A+4x0cos8-A+2x0cos8-B+2x0sin8-B/k+2x02A+2x02B 
+2cos8-tan2aA3/3+2cos8AB2— 4x0AB-2sin8AB/k-2/3+4cos8A3/3+2B3/3+2AB2+2BA2 
— 2tan2a-AB2— 2tan2a*AB2— 2cos8*tan2a/3— 2xocos8— Xotan2a— cos8*B— 2cos8*A— 2x0sin8/k+2x02 
sin8/k+2x02cos8-2x02+x0tan2a-A2+x0cos8-tan2a+x0tan2a-B2+cos8-Btan2a-2x0cos8-Btan2a+ 
3cos8-A2B— 2tan2a-A3/ 3— 2tan2a-B3/  3+cos8-B3/3+2A3/ 3-2x0B2+2x0tan2a-AB 
-2x0A2-2x0cos8-BA+x0cos8-tan2aB2-x0cos8-A2-x0cos8-tan2a-A2}  (IV-1 4a) 

This  is  the  major  result  of  our  investigation.  It  represents  the  pitch  -damping  coefficient  of  the 
lower  airfoil  in  the  presence  of  the  upper  airfoil.  Note  that  both  airfoils  are  oscillating  in  pitch. 
Also  the  two  airfoil  oscillation  may  lag  by  a  phase  angle  8. 

If  B=0  and  8—180°  the  above  result  for  the  pitch  damping  coefficient  is  the  same  as 
that  of  an  oscillating  airfoil  mounted  close  to  a  stationary  airfoil  or  wall  at  a  distance  d/2(M.  F. 
Platzer,  H.  G.  Chalkey).  As  will  be  shown  later  the  same  result  can  be  derived  using  the 
elementary  theory.  The  pitch-damping  coefficient  for  that  case  is  given  by  the  relation: 

C^-=— 2tana{x02(4— 2A)— Xo(3A2-4A+2+(A2-2)(tan2a-l))+2A3-2A— 2(2-A3)/3  (tan2a-l)] 

(IV-1 4b) 


Where  A=2hcota,  and  h  is  the  airfoil  distance. 
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For  die  special  case  of  zero  stagger  B— 0  and  out-of-phase  oscillation 
general  result  equations  (TV-12)  and  (IV-13)  simplifies  to  (IV-14). 


This  agrees  with  the  result  previously  obtained  by  Platzer  and  Chalkey  (1972)  if 
A=d.cota  in  replaced  by  A==2hcota  in  replaced  by  A— 2hcota  used  by  Platzer  and  Chalkey. 
Two  unstaggered  airfoils  oscillating  out-of-phase  are  equivalent  to  a  single  mounted  at  a 
distance  h=d/2  from  a  solid  wall  of  a  non-moving  airfoil. 

It  can  be  seen  from  Fig.  5  that  the  interference  is  highly  destabilizing  with  increasing 

amounts  of  interference. 

Note  that  A=:0.8  implies  a  larger  amount  of  interference  than  A=0.9,  for  example.  At  A—0.9 
instability  is  encountered  from  low  supersonic  flow  to  a  maximum  of  M— 3,  which  at  A— 0.8 
instability  occurs  up  to  much  larger  supersonic  Mach  numbers. 


Another  interesting  case  occurs  for  in-phase  oscillation  8=0  of  two  unstaggered  airfoils.  The 
instability  boundaries  are  shown  in  Fig.  6.  Again,  it  is  seen  that  increasing  amount  of 
interference  leads  to  increasing  regions  of  instability. 
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Fig.  6  shows  curves  of  zero  pitch  damping  coefficient  (with  torsional  stability 
boundaries)  for  an  airfoil  mounted  close  to  another  airfoil  when  B— 0  but  oscillating  in  phase. 

Following  the  same  procedure  a  special  case  of  the  above  will  be  calculated.  In  that 
case  the  upper  airfoil  is  assumed  to  be  stationary  while  the  lower  airfoil  (store  airfoil)  oscillates. 
The  results  will  be  compared  with  previous  analyses  (M.F.  Platzer,  1971). 

Starting  again  with'  Sauers  solution  for  the  slowly  oscillating  airfoil  in  unbounded 
supersonic  flow  it  is  shown  how  this  solution  can  be  extended  to  consider  the  upper  airfoil 
(wing)  interference.  To  further  simplify  the  comparison  it  is  assumed  that. 

Vm2- i— <i 

c 

Therefore,  the  moving  airfoil  is  mounted  at  half  distance  from  the  stationary  upper  airfoil  in 
compare  to  the  previous  analysis.  The  solution  is  given  by  the  following  set  of  equations. 

X(x,y)=g(z) 

vP(x,y)  =h(z)+iMyg(z)/ cosot  where  z=x— ycotoc 
VF=  h(z)+iMy  g(z)/cosa  where  z=x+ycota 

h(z)  and  g(z)  are  arbitrary  functions  for  positive  arguments  of  z  and  zero  for  z<0,  it  is  obvious 
that  the  expressions  of  the  above  functions  should  be  consistent  with  the  last  assumption. 

Using  the  two  solutions  for  left  and  right  running  Mach  waves  in  the  supersonic  flow 
field  of  slowly  oscillating  airfoils  pressure  coefficients  can  be  obtained.  In  order  to  do  so  the 
flow  field  between  the  two  airfoils  should  be  divided  into  several  zones.  The  number  of  zones 
depends  upon  A.  It  is  obvious  that  for  A>1,  there  is  one  zone  along  the  airfoil,  for  0.5<A<1 
there  are  two.  So  there  are  three  regions.  In  the  first  zone  I  there  is  no  interference  from  the 
airfoil  wing  (upper  wing),  which  means  that  the  lower  airfoil  (zone  I  only)  does  not  sense  the 
upper  airfoiL  In  that  case,  the  solution  has  the  form: 

(f>(x,y,k)  =  g(z)  +k{h(z>-  iMyg(z)/ cosa} 
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The  boundaiy  conditions  that  must  be  satisfied  at  y=0  are: 

Xy  =-cota-g'(z)=-l 

'Py=-cota-h'(z)-iMyg(z)/cosa=-i(x -Xp ) 

From  the  last  two  equations  it  is  obvious  that: 
g(z)=ztana 

h(z)=— iztana(x0+ztana/2) 

Assuming  a  slowly  oscillating  airfoil,  the  solution  to  the  PDE  has  the  form: 

<p(x,y,k)  =g(z)+k[h(z)-iMyg(z)/cosa] 

When  z  =  x— ycota 

In  order  to  calculate  the  pressure  distribution  on  the  upper  surface  of  the  lower  airfoil 
in  the  zone  I  the  following  relation  should  be  used: 

cp  =— 2[5k<p  +(pj 

Substituting  the  known  functions  (p  and  tpx 
(p=ztana  +k{-iztana(x0+ztan2a/2)-^Myztana/  cosa)} 
cpx=d/dx[ztana  +k{— iztana(x0+ztan2a/2)  — (iMyztana/  cosa)}]= 

=  tana-iktana[xo+ztan2a-(Myz/cosa)]= 

=  tana-ik[tana(x0+ztan2  a)+ (y/  sina-cosa)  (sina/  cosa)  ] = 
=tana-ik[tana(x0+ztan2a)+fy/cos2a)]= 
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=tana-ik[tana(x0  +(x— ycota)tan2a+(y/cos  a)]  — 

=tana-ik[tana(x0  +  x  tan  a-ytana+  (y/ cos  a)]  = 

=tana-ik[tana(x0+  x  tan  a-ycota] 

Therefore,  the  pressure  distribution  is  given  by  the  following  relationship. 

cp=-2{tana-ik[tana(x0+xtan2a-ycota]+ik[ztana+k{-i2tana(x0+ztan2a/2)-iMyztana/ 

cosa}]} 

Taking  into  consideration  that  higher  order  terms  of  k  are  neglected  the  pressure 
distribution  is  given  by  the  following  relationship: 

cp=2tanot{— l+ik[xo+  x  (tan2a-l)+2ycota]} 

So  at  y=0 

cp=2tana{— 1  +ik[x0+  x  (tan2a-l)]} 

In  the  zone  II: 

cp(x,y,k)  =g(z)  +k[h(z)-iM(yg(z)/cosa]+g1(zi)+k(h1(z1)+iM(ygi(z1)/cosa] 

where:  • 

Zi=x-A+ycota 

In  order  to  calculate  the  functions  g(z),h(z)  the  boundary  conditions  will  be  applied. 

The  boundary  conditions  that  must  be  satisfied  at  y— h  are. 

Xy  =-cota*g'(  z)+cota-gi'(  zx)  =  0 
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'Py=-cota-h,(z)+cota-hi,(zi)-iM[g(z)-gi(z1)]/cosa+iMAcota[g  (z)+gl  (Zl)  ]/2cosa=  0 

After  algebraic  manipulation: 

cp=-2tana{— l+ik[xo+x  (tan2a-l)]+2A-2ycota} 

For  the  region  III: 

<p(x,y,k)  =g(z) +k[h  (z)-iM(yg(z)/ cosa  ]+gi(zi)+k[hi(zi)+  iM(ygl(zi/ cosa]+ 

g2(z2)+kpi2(z2)-iMyg2(z2  )/cosa] 

where: 

z2=x— A-ycota 

The  boundary  conditions  that  must  be  satisfied  at  y-0  are: 

Xy  =-  cota-g'(  z)  +cota*gl'(  zt)  -cota*g2f  ( z$=  -1 

vPy=-cota-h'(z)+cota-h1,(z1)-cota-h2,(z2)-iM[g(z)-g1(z1)+g2(z2)]/ cosa  =  —  i(x  — xo) 

After  several  manipulations: 

cp=-2tana{3+ik[-3x0+  (3x— 2A)  (tan2a-l)]-2A/cos2a} 

Using  the  results  from  the  three  zones,  the  pressure  jump  across  the  lower  blade  can  be 
calculated  from  the  following  relationship 

Acp(x,0)=cp(x,0+)—  cp(x,0") 

In  order  to  calculate  the  cp(x,0~)  previous  results  will  be  used  (two  airfoil  oscillating  case)  for 
the  lower  surface  of  the  store  akfoil(store): 
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cp=2tana  { 1  +ik[x0+x(tan2a-l)] } 

SoAcp(x,0)=cp(x,0+)-  cp(x,Cf) 

Two  integration  intervals  can  be  distinguished  in  the  lower  airfoil  taking  into  consideration  that 
0.5<A<1: 

1)  0<x<A 

2)  A<x<2A 

In  the  first  integration  interval  is: 

AcpI(x,0)=cp(x,0+)-  cp(x,0~)= 

=-2tana{  1— ik[tana(x0+ x(tan2a— 1)} = 

=-2tana{  l-ik[tana(x0+  x  (tan2a— 1)}  and  y=0 
Acpi(x,0)=4tana{  l-ik(xo+  x  (tan2a-l)} 

The  latest  expression  is  a  complex  number  with  real  and  imaginary  part.  In  order  to  calculate 
the  pitch-damping  coefficient  only  the  imaginary  part  of  Acp  is  needed.  Therefore: 

ImAcpI(x,0)=4tana{-x0+  x  (tan2a— 1)} 

In  the  second  integration  interval: 

ImAcpn(x,0)=+2iktana[-4x0+(4x-2A)(l-tan2a)-2A/cos2a] 

The  lift  and  moment  coefficient  can  be  calculated  from: 

Q=of  *Acpdx 
C^of  1Acp(x-xa)dx 


By  definition: 
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Cjvt  0O[  C^e  +ik  Ctf] 

Q{0'=  QdPT^QaP'ii 

For  the  purpose  of  the  analysis,  only  the  pitch-damping  coefficient  is  required.  Therefore,  in 
order  to  facilitate  the  calculations  only  the  necessary  imaginary  part  of  the  will  be 
calculated. 

C„ff=  C^'j+CuQ'x  =(l/ik)JAImAcpI(x,0)  (x-x^dx+J  1ImAcpn(x,0)  (x-x^)dx= 

JA— 4tana(x0+x(tan2a— l)(x-x0)dx+Al1 2iktana[-4xo+(4x— 2A)(1— tan2a)— 2A/ cos2a]  (x— Xo)dx= 

=-2tana{x02(4-2A)-x0(3A2-4A+2+(A2-2)(tan2a-l))+2A3-2A-2(2-A3)/ 3  (tan2a-l)] 

The  above  result  is  in  agreement  with  the  previous  results  derived  by  M.  F.  Platzer  and 
H.G.  Chalkey  but  the  last  one  has  been  derived  using  die  elementary  theory. 

Therefore  the  elementary  theory  that  has  been  developed  predicts  the  above  result  as  a 
special  case  for  B=0  and  8=180°  and  represents  for  that  case  the  pitch  damping  coefficient  of 
an  oscillating  airfoil  mounted  dose  to  a  stationary  airfoil  or  wall  at  a  distance  h/2.  It  is  also  the 
same  result  of  the  analysis  of  supersonic  flow  past  a  slowly  two-dimensional  airfoil  in  a  wind 
tunnel  with  porous  wall  as  presented  by  M.F.  Platzer,  1971. 
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The  second  general  case  that  will  be  examined  in  the  thesis  is  shown  in  fig.  2.  In  that 
case,  the  Mach  cone  is  wide  enough  to  include  the  wing  airfoil 


Using  Sauer’s  solution  as  previously  and  dividing  the  flowfield  in  2ones  as  depicted  in 
Fig.7  the  perturbation  potential  <p  in  zone  I  is  given: 

<P(x,yh)  =g(z)+k[h(z)+  iMyg(z)/ cosa] 

Where: 

z=x-A-B  +ycota  at  y=h 

The  boundary  conditions  that  must  be  satisfied  at  y=h  are: 

Xy  =cota-g'  ( z)  =  -e1^ 
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.g 

*Py  =  cota-h'(  z)+iMg(z)/ cosa+iMy  g'(z)/cosa  =-iel  (x-xo-B)  : 

After  integration  and  algebraic  manipulations: 

g(z)  =-ztana*e‘° 

o 

h(z)  =iztana*elO(x0+ztan2a/2+Mh/cosa) 

Where  again  the  condition  has  been  imposed  that  the  functions  g(z)  and  h(z)  are  zero  for 
negative  and  zero  arguments. 

(|)x=-tana-e‘^+iktana-e1^  { (xo+Zitan2a+Mh/ cosa)— My/ cosa} 

=  — tanct-e*^  +iktana-e1^{(xo+(x— A— B+ycottt)tan2a+Mh/cosa)— My/cosa} 

=tana-e1^-iktana-e1^{(xo+(x— B)tan2a— A  tan2a+ycota+h/ sina-cosa) 
y/sina*cosa}  = 

=-tana-eI^+iktana*e1^  {  (xo+ (x-B)tan2a+ A— ycota} 

Because  h=A/cota ,  tana— l/cosa*sina=cota 

This  gives  the  following  result  for  the  pressure  distributions  on  the  lower  surface  of  the  upper 
airfoil  (wing  airfoil)  in  zone  I: 

cp  =-2[ik(p  +tpj 

cp=2tana-el^-2iktana-e1^  {(xo+(x-B)(tan2a-l)+2A— 2ycota}= 
cp=  2tana-e‘^  { 1— ik  {  (xq+  (x-B)  (tan2a-l)+2A-2ycota} 
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Cp=  2tana-e^{l-ik{(xo+(x-B)(tan2a-l)} 

In  order  to  calculate  the  pressure  on  the  upper  surface  of  the  wing  airfoil  the  same  approach 
■will  be  followed  and  the  result  is: 

cp—  2tana-e,^{-l+ik{(x0+(x-B)(tan2a-l)} 

To  calculate  the  Acp  of  the  upper  airfoil  (wing  airfoil)  iu  the  interval  of  zone  I  we  have  to 
calculate: 

Acp  (x^j)=:cp  (x£~)— cp  (Xjh4) 

Therefore: 

Acp=  4tana*el5  {l-ik{  (xo+  (x-B)(tan2a-l)} 

Following  the  same  arguments  as  previously  for  the  zone  II. 
z1=x— A— ycota  — B 

gi(zj)  =Zitana-(l-el5) 

The  second  boundary  condition  gives: 

-cota-h1'(2)+cota-h'(z)-iM[-g(z)+gi(zi)]/cosa=-i(x-x0-B-A) 

-cota-h1'(z)+icota-tanaei5(x0+ztan2a+Mh/cosa)-iMzitana/cosa+i(x-x0-B-A)-0 

Using  the  following  trigonometric  identity: 

» 

ztan2a-z/cos2a=-z  we  have: 
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hi(zi)  — izitana{(xo+  z1tan2a/2)(-l+  e‘^)+Mh  e^/cosa}. 

In  order  to  calculate  die  pressure  distribution  on  die  upper  surface  of  the  lower  airfoil  in  the 
zone  II  the  following  relation  should  be  used: 

cp  =— 2[ikcp  +<f>x] 

Substituting  the  known  functions  <p  and  <px 

cp=z1tana-(l-e'^)-ztana-e1^+k{iz1tana{xo+z1tan2a/2)(-l+e,^)+Mhe1^/cosa}+iztana-el^(x0+ 
ztan2a/2 +Mh/cosa)} 

<pI=d/dx[z1tana-(l-e1^)-ztana*e1^+k{iz1tana{x0+  Zj  tan2a/2)  (-1  +e'^) + Mhe^/co  sa}  +  iztana 
e1^  (x<,+  ztan2cx/2  +Mh/  cosa).}]= 

=tana-  (l-e'^)-tana-e‘^+k  (itana{  (xo+Zitan2a)  (-1  +e1^)+Mhei^/cosa}  +  itana-e1^(x0+  ztan2a 
+Mh/  cosa)} 

=  tana-(l-2e,^)+kitana[{(x0+  Zjtan^^l+e'^+Mhe^/cosaJ+e1^  (xq+  ztan2a+Mh/cosa)] 
Taking  into  consideration  that  y=0  we  have: 

=tana  {  (1— 2el^)+ki[(xo+  Z!tan2a)(— 1 +2el^)+2Mhe‘*V cosa] } 

cp  =— 2[ik(p  +cpx]=— 2[tana*(l— 2e1^)+kitana{(xo+  Zitan2a—  Zil)(— 1 +2e*^)+  2Mhe1^/cosa}] 

The  last  relation  gives  the  cp  on  the  upper  surface  of  the  lower  airfoil  (store)  in  the  zone  ILIn 
the  zone  III  a  similar  approach  will  be  used  to  calculate  the  h2(z^)  and  The  variable  z2 , 
taking  into  consideration  the  fig.  2,  is  given  by  the  following  relationship: 

Z2=x— 2A+ycota  -B 
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The  first  boundary  condition  that  must  be  satisfied  at  y=h  is: 


g2(z2)  =z2tana*(l-e‘8) 

The  second  boundary  condition  gives: 

-cota-h1'(z1)+cota-h'(z)+cota-h2'(z2)+iM[g(z)-g1(zi)+g2(za)]/cosa+iMhcota[g'(z)+gi'(zi) 

+gz  (Z2)] /cosa=-i(x-xo-B-2A)  e1 

cota-h2'(z2)+icota-tana(x0+z2tan2a)-el8Mz2tana/cosa4iMh(2-3e^/cosa4i(x-xo-B-2A)ei8 

=0 

cota-h2'(z2)=^(x0+z2tan2a)+e,^z2/cos2a-  iMh(2-3e^)/cosa-i(z2  -xo)  e"8 

h2(z2)=-iz2tana  (x0+z2tan2a/2)(l-e1^)-  iMh  tana  Z2  (2-3e^)/ cosa 
In  a  similar  manner,  as  previously,  the  cp  can  be  calculated: 

cp=z1tana*(l-et8)-ztana-e,8+z2tana-(l-el8)+k{izitana{xo+z1tan2a/2)(-l+eiS)+Mhei5/cosa} 

+i2tana*e,^(x0+ztan2a/2+Mh/cosa)+itana(x0+z2tan2ot/2)(l-e*^)-iMhtanaz2(2-3e,^)/cosa 

“iMy  fe(z)“gi  (zi) +&(z2)l  / cosa} 

(pi=tana-(2-3e‘^}+kitana[(x0+z1tan2a)(-l+e,^)+Mhe^ /cosa}  +ei8(x0+ztan2a+Mh/ cosa)+ 
(xo+  z2tan2a)(-ei8+l)-Mh(2-4ei8)/cosa] 

<pi=tana*(2-3el8)+kitana[(xo+Zitan2a)  el8-Mh(2-6ei8)/ cosa] 

In  order  to  calculate  the  pressure  distribution  on  the  upper  surface  of  the  lower  airfoil  in  the 
zone  II  the  following  relation  should  be  used: 
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cp=-2[Ik(p  +(px] 


Substituting  the  known  functions  cp  and  (px 


cp=-2iktana-[(2-3eiVz2(2-3ei5)+(xo+z2tan2a)ei5_Mh(2-6el5)/'cosa] 

=-2iktana-[(2-3ei^)+2z2(l-ei^)+(xo+(22-l)tan2a)ei^-Mh(2-6el^)/cosa] 
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V.  METHOD  OF  CHARACTERISTICS 


In  the  method  of  characteristics,  it  is  desirable  to  define  the  coordinate  system  in  such 
a  way  that  all-possible  discontinuities  could  occur  across.  Along  this  coordinate  system,  the 
equations  of  motion  of  the  flow  field  can  then  be  treated  as  ordinary  differential  equations  that 
are  solvable  by  classical  or  numerical  techniques. 

To  obtain  the  equations  in  this  coordinate  system  (the  characteristic  directions)  in  the 
(x,y)  plane,  the  equations  of  motion  are  written  in  terms  of  the  arbitrary  intersecting 
coordinates: 

(V- 1) 

And 

v (v-2) 


If  the  first  derivatives  of  the  dependent  variables,  w,v,p  with  respect  to  £  are  made 
indeterminate  across  lines  of  T|  =  constant,  and  the  first  derivatives  of  the  dependent  variables 
with  respect  to  T|  are  made  indeterminate  across  lines  of  constant,  then  any  possible 
discontinuities  in  the  first  derivatives  will  occur  across  these  lines.  These  lines  are  then  the 
characteristics  and  their  equations  are  obtained  in  the  evaluation  of  the  indetermmacies. 
Consider  first,  two-dimensional,  steady  flow.  The  governing  equations  of  motion  are  the 
continuity  equation: 

=  0  (V-3) 

dx  dy 

The  Euler  equation  becomes: 
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(V- 4) 


-du 
U - HV 

dx 


*L+L!£  =  0 

dy  p  dx 


-du  dv 

u - hv — 

dx  dy 


+I&-0 

p  dy 


(V-5) 


In  addition,  the  energy  equation  becomes: 


-  ds  8s  A 

k — + v —  =  0 
etc  dv 


(V-6) 


Along  a  streamline: 


-  dp  dp  -i-dp  dp. 

u—+v—  =  c  \u  — — - + v  — — ] 
dx  dy  dx  dy 


(V- 7) 


Substituting  the  above  into  the  continuity  equation  one  obtains: 


dw  dv 
dx  dy 


1  [-  dp  dp 

c2  L  dx  dy  J 


(V-8) 


In  terms  of  the  new  coordinates: 

pt.Zf  +Pt,VS  +V£y]pf  =-«V/>V. -^fc+V'7,k 


(V-9) 


In  similar  fashion  the  Euler  equations  become: 

pfe  +  K  +  £cP#  =  -pfc  +  K  -  Prth 

pfo  +  &P*  =  +V^K"P^ 


-PrjVx 

(V-10) 

-PrjPx 

(V-ll) 
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The  last  three  equations  form  a  system  of  three  equations  inu^,v^,p^ .  Solving  for  pq 

by  Cramer’s  rule  and  taking  into  consideration  that  is  indeterminate  across  tj  both 
determinants  in  numerator  and  denominator  should  vanish.  Therefore,  from  the  denominator 
determinant  we  have: 

fe  +v^].[(«2  -c2)£2  +2v^Jy  +  (v2  -c2)#/]  =  0  (V-12) 

The  solution  to  the  above  equation  gives  the  equations  of  all  three  characteristics  in  the 
physical  (x,y)  plane: 

—  =  — =  tan£"  (V-13) 

dx  u 

—  =  tan(f±o)  (V-14) 

dx 

Where  C,  is  the  angle  the  streamlines  makes  with  x-axis  and  cc  is  the  Mach  angle.  The  first 
equation  describes  a  streamline  while  the  last  one  describes  left  and  right  running  Mach  lines. 
The  compatibility  equation  for  u,v,p  along  ^  and  Tj  characteristics  can  be  obtained  by  setting 
the  numerator  of  p^  equal  to  zero: 


du  dv  _  1  ^  dp 

v - u — + — cota— 

dr]  dr]  p  dr] 


=  0 


(V- 15) 


The  above  relation  must  be  satisfied  along  characteristics  (£,Tj). 

According  to  Teipel  (1962)  the  unsteady  flow  over  a  flat  plate  can  be  treated  in  a  similar 
fashion.  In  order  to  facilitate  the  calculations  the  continuity  and  Euler  equations  can  be 
rewritten  using  the  local  sonic  velocity: 


(V-16) 
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So  canceling  higher  order  terms  the  continuity  equation  becomes: 


2  8c  2  dc  du  dv 
■+ - — +  c„ —  +  c_ —  =  0 


y-\8t  y- 1  "  8x  dx  dx 


(V- 17) 


In  a  similar,  fashion  the  Euler  equations  become: 


du  2  dc  du  A 

dt  y- 1  Sx  ax 


(V-18) 


av  2  ac  dv 

—  + - ca —  +  w_ —  =  0 

ar  y- 1  ax  a^ 


(V-19) 


A  system  of  three  differential  equations  has  been  formed  for  three  unknowns  u9vyc 


To  simplify  the  calculations  with  the  assumption  of  harmonic  oscillation  we  introduce  the 
amplitude  functions: 


U(x,y)eia*  =  - — — 

(V-20) 

=  -pi-— 

VM2-1  «. 

(V-21) 

yO/  \  lot  2  1  C  C 

C(x,y)e  =  — 

Vm2-1  M2  ca 

(V-22) 

Therefore,  the  continuity  and  Euler  equations  become: 

?L+M2—+ikM2C  =  0  (V-23) 

dx  dydx 


£}L+£±+ikU  =  o  (V- 24) 

dx  dx 
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dv  1  ac  ....  A 

— +  ",'■■■■  -  ■ — +  ikV  =  0 

Sx  Vm2-1  ^ 


(V-25) 


With  the  above  variables  the  compatibility  relations  for  the  above  system  of  equations  are 


fdv\ 

Kdx)s 


fdV) \  ( dC\ 


\dxj„ 


'£)  +a|V+-i- (v/2C  -£/)]  =  0  (V-26) 

\dx )£  l  M  -1  J 


./*|r__J_(M2C-C/jj  =  0  (V-27) 


rdU 

<  &  Jitr 


]  +r«v . 

J* r  V  dx  Jstr 


+  ikU  =  0 


(V-28) 


Taking  into  consideration  the  irrotationality: 


du  dv  _  q 
ay  etc 


(V-29) 


The  non-dimensional  system  of  equations  can  be  written  as: 


at/ 


+ Vm2-i— +m2 — + m2c = o 
ax  ay  ax 


(V-30) 


at/  ac  A 

- + — +ikU  =  0 

ax  ax 


(V-31) 


1=0 

ay  ax 


(V-32) 


Transforming  the  last  system  of  equations  to  the  new  coordinate  system  (^sT|)* 


f«t)  +rt_ML(c/-c)=o 

VaxJ^  vSxJ^  a/2-i 


(V-33) 
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(V-34) 


(V- 35) 


The  last  system  of  equations  can  be  solved  using  finite  differences.  The  values  of  U,V 
and  C  on  the  upper  and  lower  surfaces  of  the  airfoils  can  be  calculated  using  die  results  of  the 
Elementary  Theory  (Chapter  ID).  The  following  relations  apply: 

U  =  q>x,  V  =  tan  aq>y,  C  =  -\(px  +  ikq>\  (V-36) 

Using  the  results  of  chapter  III  the  following  results  apply  to  zone  I: 

Uj  =  tanajl  -  ik(x0  +  jctan^  a  +  >,cota)j  (V- 37) 


Vr  =  tanaf- 1  +  ik(x  -x  +  y(cota + (V- 38) 

1  [  o  cos  a  J 

Cj  =  tan  a [- 1  +  ik(x^  +  x(tzn?  a  - 1)  +  2y  cot  a)]]  (V- 39) 

The  last  three  equations  represent  the  values  of  the  parameters  on  the  upper  surface  of  the 
lower  airfoil  in  the  interval  of  the  zone  I.  The  last  equation  gives  the  pressure  distribution  and 
it  is  half  the  value  of  the  actual  pressure  distribution  coefficient.  F ollowing  a  similar  approach 
the  following  results  apply  in  the  zone  II*  for  the  pressure  coefficient: 


C  =  2{2tanaei ^  -  tana  +  ik  tana[*0  +  *(tan^  a- 1)]  - 
2 ike*$  tan  a[xQ  +  (tan^  a  - 1)(*  -  B)  +  2d  cot  or] } 


(V-40) 
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VI.  FLUTTER  ANALYSIS 


The  lift  and  moment  are  given  from  the  following  relations: 

£  =  ±p[/2}c,<ft  (VI-1) 

^  0 

* 

«  =  ipC/!}c,(jr-JT0)&  (VI-2) 

^  0 

And 


(VI-3) 


While 


(VI-4) 


(VI-5) 


And  defining: 


\cme+ikcmg]  (VI-6) 

The  last  relations  have  been  used  in  the  elementary  theory  as  shown  in  chapter  III,  and 
are  a  result  of  the  following  analysis.  To  compute  the  non  dimensional  lift  and  moment  acting 
on  the  store  airfoil  the  non-dimensional  pressure  coefficient  P  as  defined  in  problem 
formulation  last  equation  will  be  utilized  and  consequently  the  lift  and  moment  can  be  written 
as  follows: 
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0 

(VI-7) 

M  =  J  {P(*,0+ )  -  P(x,  0‘  )}>  -  x0  )dx 

0 

(VI-8) 

And: 

P(x,y)  =  yM2C(x,y) 

(VI-9) 

Where: 

P(x,y)elat  -  P~Pa>  , 

(VI-10) 

Pox 

Although  in  the  thesis  only  pitching,  oscillation  is  considered  the  more  general  case 

include  plunging  oscillation  also.  In  that  case,  Garrick  and  Rubinow’s  (1946)  method  of 

expressing  lift  and  pitching  moment  is  used: 

L  =  ~pJuJk2e"[^-(L,  +ii2)+0„<L3  +  it4)] 

^  C 

(VI-11) 

M  =  +  iM2)+ 0o(M3  +  <M4)] 

4  c 

(VI-12) 

Where: 

k  =  —  andc  =  2o 

00 

For  a  single  degree  of  freedom  as  in  that  case  the  flutter  analysis  is  much  easier,  and  the 

equation  of  motion  is: 

i,e+(i+ig)c,,0=Me 

(VI-13) 
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The  sum  of  the  moments  of  inertia  about  the  elastic  axis  is: 


Mj  =  -IQe  (VI-14) 

In  addition,  the  elastic  restoring  moment  is: 

MR=-(\  +  ig)CBe  (VI-15) 

Moreover,  the  aerodynamic  pitching  moment  about  the  elastic  axis  is: 

M,  (VI-16) 

Since: 

B  =  00eia*  (VI-18) 

6  =  -o)2d0eia*  (VI-19) 

So  the  equation  of  motion  becomes: 

em[-I9a>20o  +(l  +  ig)C06o  =  ~pJcuJk2eia,e0{M2  +iM4 )  (VI-20) 

Separating  the  equation  in  real  and  imaginary  parts  gives: 

-IgO)2  +C6  +  ^pccc4a>2M3  =  0  (VI-21) 

gCe+^pJ4eo2M4  =  0  (VI-22) 

Using  the  results  of  Garrick  and  Rubinow  (1946): 

QeX  -  pre2  +M3  =0  (VI-23) 
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M  A  +  gGlgX  =  0 


(VI-24) 


In  order  to  solve  the  single  degree  of  freedom  problem  it  is  necessary  to  obtain  values 
of  M4  as  a  function  of  k  until  the  second  equation  is  satisfied.  The  value  of  M3  for  that  k  is 
then  used  in  the  first  equation  to  determine  if  X  is  reasonable: 


X=l-^>0 

Q 


(VI-25) 


The  non-dimensional  flutter  frequency  is: 

(0F  _  1 

a>e  4x 


Accordingly  the  speed  is: 

Uf  =  1 

c(o6  kjx 


(VI-26) 


(VI-27) 


Aerodynamic  instability  (flutter)  occurs  when  M4<  0  while  the  L,  M  coefficients  are  given  by 
the  following: 


2  1  1 

L  =  L3+iL,=  — T  {jC(x, 0+)dx  -  Jc(x,0‘)tA:} 

k  o  o 


(VI-28) 


L 

M  =  M3+iM4  =-^-{|C(x,0+)(x-x0)t&-  JC(x,O“)(x-x0)<fr} 

&  n  n 


(VI-29) 


And  defining  the  following  quantity,  the  non  dimensional  pitching  moment,  the  analysis  can 
be  greatly  facilitate: 

Cm=90[cm9+ikcme.]eii>  .  (VI-30) 
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(VI-31) 


Cmff  =~\kM* 


cme  =  ~\k2Mi 


(VI-32) 
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VII.  CONCLUSIONS 


^  Itf  7;  (19?3)  “  ‘°  ^  Ae  ^offla.-ptoe 

T:  ”  ^  ^  ~  “  *  *'  -  Of  -  «*  in  L 
P  y  “  f  ^  *“*  «*  oscillate  in  pitch  wifo  m  ^ 

^;rr  *■  -**  *  - fM  ~  -  -  *  -  -  -  - 

locus  TCts°  C7  Deeded  *>  be  dlStmguished’  Qamety  “foils  with  supersonic  leading  edge 
us  su  sonic  leading  edge  locus.  In  the  former  case,  the  leading  edge  of  the  upper 

77  "  ^  ~  °f  Ma*  ~  ^  foe  lower  airfoil  ,n  foe  laner  1 

the  upper  airfoil  is  within  foe  Mach  wedge  of  foe  lower  airfoil 

pitehd  '  of  “IKCiomc  leading  edge  locus  a  general  formula  could  be  derived  for  the 

L  loX8  r  T  ^  "02  *"  «  »**  number,  pimh 

Phase  n  P  7  ^  ^  P”™01'5,  of  wo  “**•  For  the  special  case  of  out-of- 
between  foe  two  airfoih  foe  result  *tees  wifo  that  derived  by  Hataer  and 

desfohihaingefLT^y  T  “*  ^  *“  *“  ‘  «»* 

be  verified  ““  ***  <197S>  •*■*»  <*9  could 

press  TT"*  the  °f  **  SUbSOniC  ^  «*■  piobl®  *  begun  and  foe 

charaLstics  a“  ^  ^  ^  ^  <W“  “  ~  "~l  *>  ^on,  a  nmfood  of 
charactenst.es  procedure  was  described  which  could  be  used  to  compare  results. 
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